Quantum systems coupled to a structured reservoir with multiple excitations. 
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We present a method for dealing with quantum systems coupled to a structured reservoir with any 
density of modes and with more than one excitation. We apply the method to a two-level atom 
coupled to the edge of a photonic band gap and a defect mode. Results pertaining to this system, 
provide the solution to the problem of two photons in the reservoir and possible generalization is 
discussed. 



The problem of the interaction of small systems with 
structured reservoirs is of central importance to a num- 
ber of areas including nanostructurcs in semiconductors, 
atom lasers M, aspects of molecular dynamics pj and 
atoms embedded in photonic band gap materials. A fun- 
damental difficulty in the theoretical formulation of such 
problems stems from the invalidation of the Born-Markov 
approximations, essential in obtaining a master equation, 
which is the standard vehicle in the presence of smooth 
reservoirs. Models approximating some of the features 
of such reservoirs such as superpositions of Lorentzians 
can be useful, as the introduction of "pseudo-modes" can 
lead to a Markovian master equation for a system slightly 
enlarged through the introduction of the pseudo-modes 
[^-^| . Alternatively, one may introduce decorrelation ap- 
proximations in the Heisenberg equations of motion for 
the operators of interest. But in any case, one can not 
be confident of the validity and degree of accuracy of 
the approximations. For a general density of states that 
can not be modelled by a superposition of Lorentzians, 
the dynamics can only be obtained if there is at most 
one photon in the structured reservoir |7j. At this point, 
there is no generally established approach that can pro- 
vide a description of the dynamics for a general density 
of states and multiple excitations in the structured con- 
tinuum. 

It is the purpose of this paper to present such an approach 
with illustrative applications. The basic idea relies on the 
discretization of the continuum which is thus replaced in 
the formulation by a finite (but large) number of discrete 
modes. Their couplings and frequencies are chosen so as 
to model the effect of the structured continuum to the de- 
sired accuracy. The judicious choice of this parametriza- 
tion is of critical importance to the success of this idea. 
Given the discretization, the system "atom+discretized 
continuum" can be handled through differential equa- 
tions governing the evolution of the amplitudes entering 
the Schrodinger equation. These differential equations 
are then solved numerically. The discretization of con- 



tinua in other contexts pi is an established but always 
dangerous approach requiring much care as it can lead 
to unphysical artifacts. 

To introduce and demonstrate the method, we consider 
a two-level atom (TLA) coupled near-resonantly to the 
edge of a photonic band gap (PBG) P,^o|. The photonic 
band gap material has a strongly modified dispersion re- 
lation and employing the isotropic dispersion relation in- 
troduced by John and Wang , the corresponding den- 
sity of states reads 

P ( w ) = -^i=e( w - We ) (i) 

where k is a material specific constant, oj e is the band 
edge frequency and Q(x) is the Heaviside step- function. 
Clearly, the density of modes diverges at the edge fre- 
quency which invalidates the standard Born-Markov ap- 
proximations normally employed when dealing with a 
smooth reservoir. As a consequence the reservoir can 
not be eliminated. 

The idea is to replace the density of modes in Eq. (|l|) 
near the atomic transition (which for our purposes will be 
in the vicinity of the edge frequency), by a collection of 
discrete harmonic oscillators, while the rest of the mode- 
density can be treated perturbatively since it is far from 
resonance. The frequencies and the couplings of the dis- 
crete modes are chosen such that the discrete oscillators 
best model the structured continuum near the edge fre- 
quency. To this end, we write Eq. (|l|) in a differential 
form 

AN = p (lj)Alj (2) 

For AA = 1 and introducing a discrete index, we find 
Auji = l/p(u>i) and thus 

uj 1+ i = io t + Auji = LOi + l/p(u>i) (3) 

and lo\ — oj e + S where d is chosen sufficiently small 
(6 w I0~ 2 C 2 / 3 ). The coupling g r to the discrete modes 
is found by integration of Eq. (0) 
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where uj u is the upper limit of the discretized part of 
the density of states, k u is the coupling between the 
continuum mode with frequency oj and the atom and 



where C is the effective coupling 



7T ^LU — UJ e 

of the atom to the PBG structure. We thus find 



fir 



2C 



(5) 



where N is the number of discrete modes. An alternative 
form of Eq. (g) is 

= Ui-i + Auj t = cj 4 _i + 2/p{(jJi) (6) 

which is actually the form that we have used in our im- 
plementation. 

We consider a two- level atom with ground (|<?)) and ex- 
cited (|e)) states whose energy difference is hu> . The 
atom is coupled to the structured reservoir and a defect 
mode centered at the frequency w d inside the gap. The 
Hamiltonian for this system in an interaction picture ro- 
tating at the band edge frequency u> e {% — 1) and in the 
rotating wave approximation reads 

H = A a ee + A d a d a d + £ A x a\a x + g d (a d <r + + a d a ~) 

A 



+ y^,9\(a\a- + + at a ) 



A 



(7) 



where A c = w - w e , A d = uj d - uj e , A x = uj x - uj e ; 
a + = \e){g\ and er~ = \g){e\ are the atomic raising and 
lowering operators and a ee = a + a~. The field opera- 
tors (a d ,a d ) and (a\,a\) correspond to the defect mode 
and PBG reservoir, respectively, which are coupled to the 
atom via the respective coupling constants g d and g\. 
In order to demonstrate the validity of this method, we 
first present the results for spontaneous decay i.e. the 
atom is initially excited, and we neglect the defect mode 
i.e. g d = 0. Replacing the density of modes of Eq. ([!]) 
for w < w u , by a collection of discrete modes, the wave- 
function for the full system reads 
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where the amplitudes bj correspond to the discrete 
modes, while the b\ correspond to the modes with fre- 
quency lu\ > u> u , which are treated pertubatively i.e. 
they are eliminated adiabatically. 

The time evolution of the amplitudes are governed by the 
Schrodinger equation from which we obtain 

JL i 

(9) 
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a = -A o a + - 2^ 9jbj + - \^ g\b\ 
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Formal integration of Eq. (Till) gives 



b x (t) - b x (t )e l 



Ax(*-to)/< = ?± I ^' ao (i')e A A(*-*')A 



(12) 



Since these modes are strongly off-resonant i.e. A\ ^> g\ 
and for short times ao(i') remains almost constant, ao(t') 
can be replaced by do (t) . The remaining integral over the 
exponential is easily performed with the result 
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Substituting Eq. (|T^) into Eq. (||) we have 

1 1 N 2 

d = - A o a + 9jbj - E ^7 a o 

1 1 j=l A lLAX 

b 3 = 7 A i b i + T^jOo 



(13) 

(14) 
(15) 



Converting the mode sum over A into an integral from 
u> = ui u to infinity and using Eq. ([!]), we obtain 
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where for all discretized modes gj — g r as given in Eq. 
(||) . The effect of the smoothly varying part of the density 
of modes is thus to add a vacuum shift term to the equa- 
tion of motion for the upper state amplitude which effec- 
tively shifts the level down in energy and thus towards 
the band gap where it is protected from decay. This ap- 
proximation leads to a significantly reduced number of 
differential equations and the remaining amplitudes are 
distributed over a much narrower frequency interval. Be- 
yond that, the approximation also provides a surprising 
insight into the physical process, as discussed above. 
To ensure satisfactory numerical agreement with the 
known exact solution JlJ,|lJ| for this test problem, we find 
that we need at least 150 modes. In Fig. |l|, we present 
the results obtained by propagation of Eqs. ([!(]) and 
(|l7| ) . The dotted line is for a calculation with 50 discrete 
modes, the long-dashed line is for 150 discrete modes, 
and the dashed-dotted line is for 500 modes. For com- 
parison, we also plot the exact known solution [T^-[T^1 
(solid line) which shows very good agreement with the 
calculation involving 150 modes (estimated error 2%). 
The curve corresponding to 500 modes is indistinguish- 
able from the exact solution. The calculation involving 
50 modes, exhibits revivals for longer times. These are a 
consequence of the discretization; one of the dangerous 
artifacts that one must be cognizant of. Increasing the 
number of discrete modes, the revivals appear for later 
and later times. The number of modes in our calculations 
thus determines the time scale on which the propagation 
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is free of artificial oscillations, while oj u determines the 
proximity of the envelope to the correct result. This im- 
plies considerable flexibility in the method; in the sense 
that the size of the calculation can be tailored to the time 
scale, over which the behavior of the system is sought and 
the desired accuracy. 

Having demonstrated the validity of the method, we ad- 
dress now an open problem. Adding to the system de- 
scribed above a defect mode near-resonant with the atom, 
this defect mode acts as a photon source that can pump 
the atom. With one photon in the defect mode and the 
atom excited at t — 0, we have the possibility of two 
photons in the reservoir, a problem not amenable to tech- 
niques employed so far. The wavefunction for the system 
can be written 

|¥(t)) = a \e, l d , 0) + b \g, 2 d , 0) + £ b,\g, U, lj) 

3 

+ £ a, |e, 0, 1,) + hk\9, 0, l jt Ik) (18) 

3 3,k 

where the states involved are product states and, for 
instance, \g, Id, lj) = |<jr)|ld)|lj) where |l<j) is the one- 
photon state of the defect mode and is a one- 
photon state of the reservoir. The amplitudes obey the 
Schrodinger equation, and through the perturbative elim- 
ination of off-resonant modes as described above, we find 

1 / 9?N \ 
a = - A Q + A d J - a Q 

l \ U! u - LU e J 

1 1 N 

+ -V2g d b + -2_,93 h 3 ( 19 ) 

1 1 3=1 

b a = -A d b + -V2g d a (20) 

l fx x glN \ l , 

Oj- = - A + Aj 2S a 3 + -g d b 3 

! \ UJu - We / % 

1 N 1 

+ - Y. 9kb 3 k + y\ 2-y .,/.,,, (21) 

fc=l 

bj = -(Aj + A d )bj + -gja Q + -g d a 3 (22) 
bjk = r(A fe + Aj)b jk + ^-gka 3 + ^-gjak (23) 

hi = \^33 + \^!h"j (24) 

where j, k are mode indices and for all discretized modes 
9j = 9k = gr- For the purposes of this example, the 
frequency uo d of the defect mode is inside the gap as de- 
termined by the value of A^ and the atomic transition 
on resonance with uj d , i.e. A Q = A^. 
This set of equations is solved numerically with the re- 
sults presented in Fig. ^|. We plot the atomic inversion 
(solid line) , the mean- photon number in the defect mode 



(long-dashed line) and the populations in the one-photon 
sector (dot-dashed line) and two-photon sector (dotted 
line) of the reservoir Hilbert space, respectively, as func- 
tions of time. From the figure, we find that there is an 
exchange of energy (oscillation) between the defect mode 
and the one-photon sector of the reservoir. This oscilla- 
tion must involve the atom, since the defect mode is not 
directly coupled to the reservoir, but is not reflected in 
the atomic inversion. As is evident in Fig. although 
photons are exchanged between the defect mode and the 
reservoir through the atom, after some initial time, the 
atomic population remains practically constant; a rather 
surprising effect. The results presented in Fig. 2, have of 
course been tested for convergence in terms of number of 
modes, lo u , etc. 

As the defect mode is pushed further into the gap, we 
find that the oscillations of the atomic population begin 
to extend to increasingly longer times. On the contrary, 
a change in the magnitude of g d in relation to C 2 / 3 does 
not seem to affect the atomic oscillations for longer times, 
but it does affect the relative oscillations of the excita- 
tions in the defect mode and the reservoir, as we will 
discuss in detail elsewhere. 

In conclusion, we have developed an approach that is ca- 
pable of providing solutions to a class of problems which 
can only be treated approximately through other tech- 
niques. It is applicable to small systems coupled to a 
density of modes of any form and has allowed us to solve 
problems involving multiple excitations in the contin- 
uum. In addition to the implementation outlined here, 
we have explored various other forms of discretizations, 
as well as other densities of modes, with good agreement 
with other exact results in those cases that are available. 
This demonstrates the generality and versatility of the 
approach, which could be readily employed in other con- 
texts such as, for example, waveguides where the density 
of modes is also singular. Jl^Jl7[] 

We demonstrated, in addition, that off-resonant contin- 
uum modes can be eliminated perturbatively. The effect 
of this approximation is to reduce the number of differen- 
tial equations to be propagated, leading thus to a drastic 
enhancement of the computation speed, essentially with- 
out compromise in accuracy. The number of equations 
to propagate scales roughly as N p , where N is the num- 
ber of discrete modes and p the number of excitations. 
The ultimate limitation of the method is determined by 
computer memory as demanded by each problem. In par- 
ticular, the study of multiple excitations will probably for 
the time being be limited to 4 or 5. For the purpose of 
presenting the method, we limited our discussion here to 
two photons. Results from work for more photons will 
be presented elsewhere, but we do mention here that we 
have also obtained fully converged results for three pho- 
tons in the reservoir. 

Providing solutions for the dynamics of the system is one 



3 



aspect of this approach. Perhaps an equally important 
aspect is the insight gained by the possibility to combine 
the perturbative treatment with the non-perturbative 
treatment of the rest. As discussed above, this sheds light 
on the physical effect of the modes around the band edge 
as compared to the smooth distant part of the density 
of modes. In addition, the possibility to monitor the dy- 
namics of additional photons in the structured reservoir 
may prove very valuable when considering the validity of 
approximations necessary in other schemes. 
It should be mentioned in closing that a recently pro- 
posed formal approach to similar non-Markovian prob- 
lems, based on the quantum state diffusion formalism, 
has been presented by Diosi, Strunz and Gisin |Q. At 
this point we are aware of the application of the method 
to a relatively tractable problem involving a standard 
cavity reservoir. Its potential, however, does not seem at 
first sight to be limited and it will be interesting to see 
and explore its applications to problems involving more 
complicated densities of modes. 
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FIG. 1. The population in the excited state as function of 
time. The solid line is the exact solution. The dotted line 
is for N = 50. The long-dashed line is for N = 150 and 
the dotted-dashed line is for N = 500. The insert shows a 
close-up of the long-time behavior. Parameters: A D = and 
9d = 
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FIG. 2. The evolution of the system is plotted as func- 
tion of time. The solid line is the population in the upper 
atomic state. The long-dashed line is the mean photon num- 
ber in the defect mode, the dot-dashed line is the population 
in the one-photon sector of the reservoir Hilbert space and 
the dotted curve is the population in the two-photon sec- 
tor of the reservoir Hilbert space. Parameters: TV = 150, 
g d = C 2 / 3 , A„ = A d = -0.1C 2 / 3 
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